00 

o 

O 
(N 

D 
00 



D 
u 

-)— » 

a 

a 
o 

(N 
> 

00 

oo 

p 
o 

oo 

o 



% 



Vector chiral and multipolar orders in the spin-1/2 frustrated ferromagnetic chain 

in magnetic field 
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We study the one-dimensional spin-1/2 Heisenberg chain with competing ferromagnetic nearest- 
neighbor Ji and antiferromagnetic next-nearest-neighbor J2 exchange couplings in the presence 
of magnetic field. We use both numerical approaches (the density matrix renormahzation group 
method and exact diagonahzation) and effective field-theory approach, and obtain the ground-state 
phase diagram for wide parameter range of the coupHng ratio J1/J2. The phase diagram is rich and 
has a variety of phases, including the vector chiral phase, the nematic phase, and other multipolar 
phases. In the vector chiral phase, which appears in relatively weak magnetic field, the ground state 
exhibits long-range order (LRO) of vector chirality which spontaneously breaks a parity symmetry. 
The nematic phase shows a quasi-LRO of antiferro-nematic spin correlation, and arises as a result 
of formation of two-magnon bound states in high magnetic fields. Similarly, the higher multipolar 
phases, such as triatic [p — 3) and quartic {p — 4) phases, are formed through binding of p 
magnons near the saturation fields, showing quasi-LRO of antiferro-multipolar spin correlations. 
The multipolar phases cross over to spin density wave phases as the magnetic field is decreased, 
before encountering a phase transition to the vector chiral phase at a lower field. The implications 
of our results to quasi-one-dimensional frustrated magnets (e.g., LiCuV04) are discussed. 

PACS numbers: 75.10.Jm, 75.10.Pq, 75.40.Cx 



I. INTRODUCTION 

There is resurgence of theoretical interest in the one- 
dimensional frustrated ferromagnetic Heisenberg model 
in magnetic field. 
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where the nearest-neighbor exchange is ferromagnetic 
Ji < 0, the competing next-nearest-neighbor exchange 
is antiferromagnetic J2 > 0, and s; is a spin-^ oper- 
ator on the site I. The model has recently attracted 
much attention as it is considered to describe magnetic 
properties of quasi-one-dimensional edge-sharing chain 
cuprates, such as Rb2Cu2Mo30i2 (ReLiJu) and LiCuV04 
(Ref. 111) . In particular, there have been intensive exper- 
imental studies of LiCuV04, exploring an unusual phase 
transition in magnetic field from a spiral-ordered phase 
to a modulated-coUincar-ordered phasc^^-^'^ and a multi- 
ferroic behavior J^ii^ii^ 

From a theoretical point of view, the J1-J2 spin chain 
(H]) is of special interest as it is the simplest of the frus- 
trated quantum spin models and provides a good testing 
ground to look for exotic quantum phases induced by 
frustration. The theoretical studies over the past several 
decades have mostly considered the case where both cou- 
plings are antiferromagnetic, Ji > and J2 > 0. It has 
been established that in zero magnetic field the ground 
state of the antiferromagnetic J1-J2 spin chain under- 
goes a phase transition from a critical phase with gap- 
less excitations for J2 < J2c = 0.2411Ji to a gapped 
phase with spontaneous dimerization for J2 > J2C as J2 



increases! "'^''''^^i^^i^^i^^'^^i^^ It has also been revealed that 
the model exhibits cusp singularities and a 1/3-plateau 
in the magnetization curv o^^'^^ as well as a vector chiral 
order in the case of anisotropic exchange couplings^^i^li^ 
or under magnetic £ieldj ^^'^°i'^^'^^ 

In this paper we concentrate on the ferromagnetic 
case ( Ji < 0) of the J1-J2 spin chain ([1]) in magnetic 
field which partially polarizes spins to the -\-z direction. 
We show that the ground-state phase diagram in the 
case is a zoo of exotic quantum phases, using the nu- 
merical density-matrix renormahzation group (DMRG) 
method^^i^ii^i^ exact-diagonalization method, and ef- 
fective field theories. We find a phase with long-range 
vector chiral order and phases with various kinds of mul- 
tipolar spin correlations, most of which have not been 
known to appear in this model. 

Let us briefly review established results from previous 
studies on the ferromagnetic J1-J2 spin chain and intro- 
duce our new findings. 

In zero field the ground state is ferromagnetic for 
Ji/ J2 < —4 and spin singlet for —4 < J1/J2 < 0; the 
nature of the spin singlet ground state is not well under- 
stood. The ground state manifold has extensive degener- 
acy at the phase boundary J1/J2 — — 4i^i^ In magnetic 
field the spins order with a helical magnetic structure 

si /s — (sin 9^ cos 0^ , sin 9^ sin 0^ , cos 9^) (2) 

in the classical limit (s = |s| 3> 1), with a pitch angle 

(/,'= = (1)1^^ - 0? = ± arccos(- J1/4J2) (3) 

and a canting angle 



9" = arccos[4/iJ2/s(Ji + 4J2)^], 



(4) 



when —4 < J1/J2 < 0. One might expect that this he- 
Hcal magnetic order should be completely destroyed by 
quantum fluctuations in the quantum limit s = 1/2. It is 
important to notice, however, that a part of the broken 
symmetries in the classical helical spin configuration may 
remain to be spontaneously broken even in the quantum 
limit. Indeed, the chirality (the sign of 0^) of the heli- 
cal spin configuration is Z2-valued and can be broken in 
(1+1) dimensions. The chirality can be measured with 
the vector chiral order parameter 



iO 
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(5) 



with ', 

9 • 2 
s sm 



— 1 and 2; its classical value is 



,(") 



sin(n(/)'^). In this paper we show, for the first 

(n) 

time, that the vector chirality kI is long-range ordered 
in the weak-field regime of the phase diagram of the fer- 
romagnetic J1-J2 model^S. We also show that the vector 
chiral order parameters satisfy the relation 



Ji(k 



(i)\ 



2J2(Kp^) = 0, 



(6) 



where (■ ■ •) denotes average in the ground state. This 
implies that the spin current, Jij{si x SjY, flowing on 
the link connecting the sites i and j {Jij = Ji or J2) 
is confined and circulating in each triangle made of the 
three neighboring sites. Incidentally, we note that the 
classical helical configuration ^ satisfies Eq. ([6]). 

The vector chirality ([5]) is an antisymmetric prod- 
uct of two spin-i operators. This is an example of 
the p-type nematic operator introduced by Andreev and 
Grishchuk.'*'^ In this paper, we shall reserve the term 
"nematic" for symmetric products (termed n-type in 
Ref. \4u) and call the antisymmetric product ^ the vec- 
tor chirality. Examples of what we call nematic operators 



^x^—y^ '^i ^7 



s^4' 



ixy 



= <4 






(7) 



which can be thought of as members of quadrupolar spin 
operators. 

Interestingly enough, the phase diagram of the fer- 
romagnetic J1-J2 spin chain has a Tomonaga-Luttinger 
(TL) liquid phase with quasi-long-range antiferro- 
nematic order Q = Q^2_,i2 — i 



ix'^ —y^ 



s, s,- , where 



s~ = s^ — is^,: As first pointed out by Chubukovji this 
nematic order is realized due to pairing of two magnon ex- 
citations. The paired magnons are the low-energy excita- 
tions of the TL liquid with the nematic quasi-long-range 
order. This was confirmed recently by numerical calcu- 
lation of nematic correlation function at J1/J2 = — li^ 
In this paper we explore wider region of the parameter 
space and show that the nematic TL liquid phase occu- 
pies a large part of the phase diagram. 

One can generalize the quadrupolar spin orders to 
higher multipolar orders. For example, one can de- 
fine ocutupolar triatic order"*^ O = O^^^^xy^ + 



y'^~3x^y 



St Sj Sj 
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and, similarly, the hexadecapolar order iJ_ 



i/.4_ 



x^—Gx'^y'^+y^- 



- iH,r 



Sj s Sj, S; , which we dub 



the "quartic" order, and so on. In fact, signatures of 
binding of three or four magnons are found in recent nu- 
merical studies of magnetization curvea^i^ and of mul- 
timagnon instabilities at a saturation field.— In this pa- 
per we establish the existence of TL liquid phases with 
the triatic and quartic orders through the DMRG calcu- 
lation of correlation functions. It is interesting to note 
that quasi-long-range molecular superfluid phases (called 
trionic and quartetting phases), similar to the above- 
mentioned triatic and quartic phases, have recently been 
found in a model of one-dimensional multicomponent 
fermionic cold atomsi^i^i^ 

This paper is organized as follows. In Sec. |TT] we 
present the ground-state phase diagram and briefly de- 
scribe properties of the phases newly identified in the 
present work. These are vector chiral, nematic, incom- 
mensurate nematic, triatic, quartic phases, and spin den- 
sity wave phases (SDW2 and SDW3). This section gives a 
summary of the main results of the paper. In Sec. IIIII we 
consider formation of multimagnon bound states which 
destabilizes the fully polarized state. The results of this 
consideration allow us to determine phases emerging just 
below the saturation field. In Sec. lIVI we study magneti- 
zation process of the model ([I} for several values of the 
ratio J1/J2, and find a transition from a single-spin flip 
process to a multispin flip process. The transition point 
is identified as the boundary of the vector chiral phase 
in the phase diagram. The remaining sections present 
detailed analysis of correlation functions in each phase, 
which we calculate using the DMRG method. In Sec. |Vl 
we consider the vector chiral phase. After a brief review 
of bosonization approach due to Kolezhuk and Vekua^ 
which is valid for |Ji| <C J2, numerical results of the 
DMRG calculation are presented. In Sec.|Vl]we discuss 
the nematic phase. We introduce a hard-core bose gas 
of two-magnon bound states as an effective theory for 
the nematic phase. We find good agreement between the 
theory and numerics of various correlation functions in 
the nematic phase. In Sec. IVIII we show our numerical 
results for the incommensurate nematic phase, which ex- 
hibits quasi-long-range order of the nematic correlation 
with an incommensurate wave number. In Sec. IVIIII we 
apply the hard-core boson theory to the triatic and quar- 
tic phases. We show that these phases can be understood 
as TL liquids of hard-core bosons which correspond to 
three- and four-magnon bound states, respectively, just 
as the nematic phase is a TL liquid of two-magnon bound 
states. We conclude with some remarks in Sec. IIXI The 
relation ^ is derived in the Appendix. 
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FIG. 1: (Color online) Magnetic phase diagram of the spin- 
1/2 zigzag chain with ferromagnetic Ji and antiferromagnetic 
J2 (a) in the J1/J2 versus h/J2 plane and (b) in the J1/J2 
versus Af plane. Crosses show the transition and crossover 
points obtained from the magnetization curves and corre- 
lation functions. In (a), symbols "VC", "N", "IN", "T", 
"Q", and "F" indicate the vector chiral (ASfot ~ 1), ne- 
matic (ASfot ~ 2), incommensurate nematic (ASfot = 2), 
triatic {ASt^t = 3), quartic (ASfot = 4), and ferromagnetic 
(fully polarized) phases, respectively. Here AS'fot is the unit 
of changes in the total Stot = X^j sf when the magnetic field 
h is swept. There are also two kinds of spin-density wave 
phases: SDW2 (A^t'ot = 2) and SDW3 (ASfot = 3), which 
are related to the nematic and triatic phases, respectively. 
The solid curve shows the saturation field /is and dotted lines 
are the guide for the eye. In (b), symbols indicate parame- 
ter points for which their ground-state phase is identified by 
analysis of correlation functions. Shaded regions in (b) corre- 
spond to the magnetization jump at the first-order transition; 
see Sec. HV] 



II. PHASE DIAGRAM 

The ground-state phase diagram obtained in the 
present work is summarized in Fig. [1] in the planes of 
(a) J1/J2 versus /1/J2 and (b) J1/J2 versus the mag- 
netization per site M. The phase diagram contains (at 
least) eight phases: vector chiral (VC) phase, nematic 
(N) phase, incommensurate nematic (IN) phase, triatic 



(T) phase, quartic (Q) phase, two kinds of spin density 
wave phases (SDW2 and SDW3), and ferromagnetic (F) 
phase. Brief explanation of these phases is given below. 
More detailed discussions on each phase will be given in 
Sees. rvHVIIIl Figure [2] shows typical spatial dependence 
of various correlation functions in these phases. 

Ferromagnetic phase: In the ferromagnetic phase, 
spins are fully polarized, M = 1/2. This phase is sta- 
ble when J1/J2 < —4 or when large enough magnetic 
field is applied. We will discuss in Sec. IIIII magnetic in- 
stabilities along the phase boundary of the ferromagnetic 
phase. 

Vector chiral phase: The vector chiral phase appears 
in small magnetic field. This phase is characterized by 
long-range order of the vector chiral correlation ([5]) . The 
ground state breaks a Z2 symmetry, as Eq. ([5]) indicates 
that the parity about a bond center is broken sponta- 
neously. We can also regard this Z2 symmetry breaking 
as choosing one of the two possible directions of circu- 
lation of spontaneous s^-spin current flow. A schematic 
picture of the vector chiral order and circulating spin cur- 
rent in the Z2-symmetry broken state is shown in Fig. [31 
where the spin chain is drawn as a zigzag ladder. Numer- 
ical evidence for the long-range order will be presented 
in Sec. |Vl Another important feature of the vector chi- 
ral phase is that the transverse spin correlation {sQs'f) is 
incommensurate with the lattice and stronger than the 
longitudinal correlation (sgsf) — {sQ){sf). 

Nematic/ SDW2 phases: At higher magnetic field up 
to the saturation field, the nematic/SDW2 phasesi'^'i'^ 
appear at J1/J2 > —2.7. These phases are a TL liq- 
uid of hard-core bosons which are actually two-magnon 
bound states with total momentum k — tt. The boson 
creation operator bj corresponds to s/sj,^, and the bo- 
son density n; = 6j6; oc ^ — sf. Since breaking a two- 
magnon bound state costs a finite binding energy, the 
transverse spin correlation {sQsf) is short-ranged, where 
■^d^ = '^0 + *'^o- Being a TL liquid, the ground state ex- 
hibits power-law decaying correlations of the single-boson 



propagator, (69 6;) 



ex (s^s^si 



S;_|_j^), and the density fiuc- 



tuations, (noui) — (uq) (ni) ex (so^f) — (so)(sf ). When the 
boson propagator decays slower than the density-density 
correlation, it is appropriate to call this phase the (spin) 
nematic phase. In the opposite case when the latter in- 
commensurate density correlation is dominant, we call 
this phase the spin density wave (SDW2) phase. The 
SDW2 phase is extended to the antiferromagnetic side 
Ji > across the decoupled-chain limit Ji — 0; it is 
called even-odd phase in Ref. [23. The boundary between 
the SDW2 phase and the nematic phase is shown by a 
dotted line in Fig. [T] 

In the semiclassical picture we can write s/ — e"*"*"' , 
where (jji is the angle of the two-dimensional vector 
(sf , s^) measured from the positive x direction, < 0/ < 
27r. The product sr^z+i = e^'^'^'+'^'+i' can be rep- 
resented by the vector Ni^i — (cos$;^2 7sin<i>i^2) with 
$i^2 = —{(f>i + (f>i+i)/'^- We now realize that we need to 
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(c) J/J2 = -2.0, M = 0.2 (SDW2) 
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(f) J1/J2 = -3.0, M = 0.2 (SDW3) 



FIG. 2: Typical behaviors of various correlation functions in (a) the vector chiral (VC) phase, (b) the nematic (N) phase, (c) 
the SDW2 phase, (d) the incommensurate nematic (IN) phase, (e) the triatic (T) phase, and (f) the SDW3 phase. Absolute 
values of spatially averaged correlation functions are plotted. (For the averaging procedure, see Sec. lVBl ) In (b) and (d), the 
triatic correlation function {s^s'^+iSJ.2^'i''^7'+i^7'+2) ^^ smaller than 10~^. 




FIG. 3: Schematic picture of the vector chiral order. The 
arrows on bonds indicate breaking of the parity symmetry by 
the vector chiral order Kj — {si x s;+„)^, whose expectation 
values obeys the relation .Ji{k^^') +2J2{k, ) ~ 0. The circu- 
lation of the s^ spin current, shown by the dashed arrows, is 
alternating, and there is no net spin current flow. 



identify Ni^i with — 7V;^ 1 because of the physical iden- 
tification ((/)/, 0/+i) = {4>i + 27r, (^/+i) = {(j>i,(j>i+i + 27r). 
We can thus consider Ni^i as a director representing the 
nematic order. We will show in Sec. IVII that the nematic 
phase has antiferro-nematic quasi- long-range order of the 
director, as shown schematically in Fig. |31 The ground 
state is not dimerized in this phase, as opposed to the 
initial proposal of Chubukovji' 

Incommensurate nematic phase: The incommensurate 
nematic phase occupies a very small region in the phase 
diagram. This phase has quasi-long-range order of the 
nematic correlation with an incommensurate wave num- 
ber. The correlation is due to two-magnon bound states 
with momentum k — tt + 6 and tt — S, instead of fc = tt in 
the nematic phase. Schematic pictures of the incommen- 
surate nematic order are depicted in Fig. [51 where the 
upper and lower pictures represent the nematic correla- 
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FIG. 4: Schematic picture of antiferro-nematic quasi-long- 
range order in the nematic phase. Ellipses represent directors 
of the nematic order on each bond. 



tion with wave number k — n + 6 and tt — (5, respectively. 
If the densities of paired magnons with k = ir + S and 
TT — 5 are different, one of the two correlation patterns in 
Fig. E] becomes dominant, and the Z2 chiral symmetry is 
broken spontaneously, as suggested by Chubukovii How- 
ever, we have found no signature of long-range order of 
the chiral correlation in our numerical calculation, which 
we discuss in Sec. IVIII 

Triatic and SDW3 phases: The triatic phase exists 
below the saturation field and next to the incommen- 
surate nematic phase. The triatic/SDWa phases are 
a TL liquid of bosons which represent three-magnon 
bound states with total momentum k = n. In anal- 
ogy with the nematic phase, the triatic order has a 
simple semiclassical picture. Writing the bound three 



magnons as s 



-i{<l>i+<pi + i+rpi+2) — g3i*i,3 yifQ 



may consider the triatic order as ordering of the angle 
^1,3 = -(</>; + 4>i+i + <t>i+2)l'i, which has the property 
$(^3 = $;.3 + 27r/3. A schematic picture of the triatic or- 
dered state is shown in Fig. [Sj In the triatic phase, cor- 
relation functions probing three-magnon bound states. 






FIG. 5; Schematic pictures of incommensurate nematic 
quasi-long-range order in the incommensurate chiral nematic 
phase. Ellipses represent directors of the nematic order on 
each bond. The numerical results in Sec. IVIII indicate that 
this chiral symmetry is not broken in the incommensurate 
nematic phase; the ground state is given by equal superposi- 
tion of the upper and lower configurations. 



TABLE I; Number of magnons, p, and total momentum k of 
the multimagnon bound states which become gapless at the 
saturation field. 



parameter range 


P 


k 


-2.669 < J1/J2 < 


2 


TT 


-2.720 < J1/J2 < -2.669 


2 


n±5 {S>0) 


-3.514 < J1/J2 < -2.720 


3 


TT 


-3.764 < J1/J2 < -3.514 


4 


TT 


-3.888 < J1/J2 < -3.764 


5 


n 


-3.917 < J1/J2 < -3.888 


6 


n 


-4 < J1/J2 < -3.917 
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-K 




FIG. 6: Schematic picture of antiferro-triatic quasi-long- 
range order in the triatic phase. Solid triangles represent spin 
structure of the triatic order formed by three s = 1/2 spins 
on each plaquette. 



such as {sq Si $2 s'[ s'i^j^-^s'^^2) ^ exhibit quasi-long-range or- 
der (power-law decay). In contrast, both the trans- 
verse spin correlation {s^s^) and the nematic correla- 



tion (ststs, s 



are short-ranged, because of a finite 
energy cost for breaking a three-magnon bound state 



'=1 ''; ^z+i/ 



in our previous studyi^ 

Inside the ferromagnetic phase in magnetic field, there 
is a finite energy gap between the ground state and ex- 
cited states. With decreasing the magnetic field, the gap 
becomes smaller and eventually vanishes at the boundary 
of the ferromagnetic phase. We define the saturation field 
hs as the magnetic field h at which a branch of excitations 
first becomes gapless as the field h is reduced. In the fer- 
romagnetic J1-J2 spin chain, the excitation mode that 
becomes gapless (soft) at h — hs is a multimagnon bound 
state. Below the saturation field the soft multimagnon 
bound states proliferate. As a result the ground state 
can change into a TL liquid with the correlation that is 
represented by the soft, bound multimagnon mode. It is 
therefore important to find out which branch of multi- 
magnon bound states is the soft mode. 

We calculate energy of p-magnon excitations using the 
method we introduced in Rcf. 6. The number of magnons 
p and the total momentum k are good quantum numbers 
of the Hamiltonian ([l]). We thus expand eigenstates in 
the sector of p magnons with the basis 



SI p 



The longitudinal spin correlation (snsf) — (sn)(sf ) shows Ink-ir, r ,i\ — V^ TT 



algebraic decay, as sf is related to the boson density. 
When the most slowly decaying correlation is the lon- 
gitudinal spin correlation, we call the phase the SDW3 
phase. The boundary between the triatic phase and the 
SDW3 phase is shown by a dotted curve in Fig. [TJ The 
detailed discussion of the correlation functions will be 
given in Sec. IVIIII 

Quartic phase: The quartic phase is a TL liquid phase 
of four-magnon bound states with momentum k = n. 
Its properties can be easily deduced by straightforward 
generalization from the triatic phase. 



III. MULTIMAGNON INSTABILITY 

We begin our study of the phase diagram by examining 
instabilities of the fully polarized state. To that end, we 
numerically calculate energy dispersion of low-energy ex- 
citations with a small number of magnons (down spins). 
The analysis presented here extends the result reported 



Vn 



e^kU/Pg- 



= 1 n=l 



where 



hi ^ 



n-l 



|FM), (9) 



(10) 



Here |FM) is the fully polarized state (s+|FM) ^ 0), n 
the system size taken to be $7 — > 00, and r^ (1 < i < p— 1) 
is the distance between the i-th and (i + l)-th magnons. 
The periodic boundary condition is imposed in this cal- 
culation. We take r.; to be in the range of 1 < r^ < rmax, 
where r,„ax is chosen so that wave function vectors can 
be stored in the computer memory. (The finite value 
of Tmax hmits the accuracy of energy calculations. For 
tightly bound magnons, errors caused by this approxi- 
mation can be made negligibly small, on the order of 
exponentially decaying tails of their wave function.) 

We numerically diagonalize the Hamiltonian matrix 
expressed in the restricted Hilbert space (r.^ < rmax) and 



obtain the lowest energy as a function of the total mo- 
mentum k for each p magnon sector. In this way we ob- 
tain energy dispersion of p-magnon bound states. In our 
previous studyS- we calculated energy dispersion of multi- 
magnon excitations for up to p = 4. Here we extend the 
calculation to include more magnons (pmax = 8), taking 
the maximum distance r,„ax to be at least 42/(p— 1). We 
thereby identify soft multimagnon modes and determine 
the saturation field hs at each value of the ratio J1/J2 
(-4 < J1/J2 < 0). 

Table U summarizes the number of magnons, p, and the 
momentum k of the multimagnon modes {p < 8) which 
become gapless at the saturation field. We note that 
the gapless modes with p < 4 in Table [J are soft modes 
giving rise to multipolar TL liquids, since non- negative 
excitation energy of 2p-magnon modes indicates a repul- 
sive interaction between bound p-magnons (see the dis- 
cussion at the end of this section). We thus find that, 
as Chubukov first pointed out,^ the two-magnon bound 
state with fc = tt is the soft mode when —2.669 < J1/J2 < 
0. Its exact wave function can be easily obtained; it turns 
out that the bound-state wave function at A; = tt has am- 
plitudes only for odd integer values of ri, which means 
that the magnons forming a bound pair are on different 
legs of the zigzag ladder. The soft mode signals emer- 
gence of a nematic phase below the saturation field. 

In the narrow range —2.720 < J1/J2 < —2.669, the 
soft two-magnon bound state has an incommensurate 
momentum fc 7^ tt. Our numerical estimate of the 
commensurate-incommensurate transition point is con- 
sistent with the exact result, {Ji/J2)c = —2.66908....- 
Beyond the commensurate-incommensurate transition 
point, the total momentum fc changes continuously as 
fc/TT = 1 - 0.67^(Ji/J2)c- J1/J2; see Fig. [H This sug- 
gests continuous nature of the transition between the 
commensurate and incommensurate nematic phases at 
h — hs- 

As the ratio J1/J2 is changed towards the end point at 
J1/J2 — —4, the magnon number p of the lowest bound- 
magnon branch increases; see Table I. The total momen- 
tum of the bound-magnon mode is always at fc = tt, 
except in the narrow region of the incommensurate two- 
magnon bound states mentioned above. We expect that 
many- magnon bound states, formed by more than seven 
magnons, should appear as J1/J2 -^ —4. In our numeri- 
cal calculation, eight-magnon bound states did not come 
down as the lowest state, which we suspect was due to 
finite-size effects coming from small rmax- 

To demonstrate stability of the multimagnon bound 
states, we show in Fig. |S] dispersion curves of bound- 
magnon excitations, as well as lower edges of continuous 
spectra of magnon scattering states, at the saturation 
field h — hs for Ji/ J2 = —3.0 and —3.6. The p-magnon 
scattering states are constructed from a set of one-, two-, 
..., {p— l)-magnon (bound) states, in total oip magnons. 
At J1/J2 — —3.0 [Fig. [HIa)], the three-magnon bound 
state is gapless at fc = tt. The branches of one, two, 
four, and five magnons have finite excitation gaps. This 
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FIG. 7: Dependence on J1/J2 of the center-of-mass mo- 
mentum k for the two-magnon bound state. The momen- 
tum deviates continuously from tt at J1/J2 — —2.669. The 
incommensurate momentum k is fitted well to k/n = 1 — 
0.67-1^—2.669 — J1/J2 as shown by the solid curve. The 
dashed line is the classical estimate k = 2arccos ( — J1/4J2) 
for two scattering magnons. 



feature is consistent with a finite binding energy of the 
three-magnon bound state. Furthermore, the state with 
the lowest energy (at fc = 0) in the six-magnon sector 
belongs to the continuum of scattering states formed by 
a pair of three-magnon bound states. This indicates that 
three-magnon bound states are interacting repulsively 
with each other. The repulsive interaction rules out the 
possibility of a metamagnetic transition (magnetization 
jump) at the saturation field, and instead induces a con- 
tinuous transition to the triatic phase which we discuss 
in more detail in Sec. lVIlIl The multimagnon dispersions 
for J1/J2 = —3.6, shown in Fig. [HJb), can also be under- 
stood in the same fashion. Here it is the four-magnon 
bound states that become gapless at the saturation field. 
The instability of the fully polarized state is driven by the 
four-magnon bound states with mutual repulsive interac- 
tions, which condense to form a TL liquid with quartic 
order, as we show in Sec. IVIIIl 



IV. MAGNETIZATION CURVE 

Having identified the soft modes at the saturation field 
h = hg, we now study magnetization process of spin 
chains of finite length, which we obtain numerically for 
various values of the coupling ratio Ji/ J2. The numerical 
results help us to deduce overall structure of the magnetic 
phase diagram. Previous studies^i^ have found that, near 
the saturation field, the total magnetization iSf^t — J^i ^i 
changes in units of AS^^^ =2,3, and 4 for J1/J2 > —2.6, 
—3 < J1/J2 < —2.8, and J1/J2 = —3.75, respectively. 
The multispin flip is a natural consequence of the forma- 
tion of stable multimagnon bound states. At lower fields 
S'tot changes by AS'^^^^ = 1. 

We obtain magnetization curves from the following 
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FIG. 8: (Color online) Dispersion curves of multimagnon 
bands at the saturation field for (a) J1/J2 = —3.0 and (b) 
Ji/ J2 ~ —3.6. The solid curves are the dispersions of bound 
states and the dashed curves show the lower edges of the 
continuum of scattering states. For clarity, only the states of 
up to 6-niagnons are shown in (a) and up to S-magnons in (b). 
The numbers printed beside the curves denote the number of 
magnons. Bound states inside the scattering continuum are 
not shown. 



procedure. With the DMRG method we calculate the 
lowest energy Eo{M) of the model ([1]) at ft, = in each 
Hilbert subspace of magnetization per site M = S^^^/L, 
where L is the number of total sites. The magnetization 
curve M{h) is then obtained by finding the magnetiza- 
tion M which minimizes Eq{M) — hM for given h. We 
have performed the calculation for open chains of up to 
L = 168 sites. We kept up to 350 states in our DMRG 
calculation. 

Figure [HI shows representative magnetization curves 
calculated at various values of Ji/J2- We clearly see that 
the total magnetization S^^^ changes by AS^q^ = 1 at 
low magnetic fields, while it shows multispin flip process 
AS't^ot > 2 at higher fields. The magnetization change 
in the high-field regime is AS'^^t — 2 for J1/J2 > —2.7, 
A5t^„t = 3 for -3.4 < J1/J2 < -2.75, and /:^S'{^^ = 4 
at J1/J2 = —3.6. These features of the magnetization 
process are consistent with our finding of stable multi- 
magnon bound states discussed in Sec. IIIII The critical 
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FIG. 9: Magnetization curves for (a) J1/J2 = —2.0, (b) 
J1/J2 = -2.4, (c) J1/J2 = -2.5, (d) J1/J2 = -3.0, (e) 
J1/J2 = —3.4, and (f) J1/J2 = —3.6. The dotted lines repre- 
sent the boundaries of the regions of AS'tot = 1 and AStJot > 2. 



field he and the critical magnetization Mc, at which the 
magnetization step changes from AS"^ = 1 to /S.S^ > 1, 
are plotted in the phase diagram shown in Fig. [TJ 

The magnetic phase diagram (Fig. [T]) has four dis- 
tinct regions characterized by AS^^^ = 1,2,3, and 4. 
These regions correspond to the vector chiral, nematic 
or SDW2, triatic or SDW3, and quartic phases, respec- 
tively. We will discuss each region in detail in the follow- 
ing sections. The phase boundary between the region of 
AS^^^ = 1 and that of AS^^^ — 2 begins from the critical 
coupling J1/J2 = —2.72 at h — hg, which is the phase 
boundary between the incommensurate nematic and tri- 
atic phases,— and appears to go towards small IJ1I/J2 
region as h is decreased. 

When —2.72 < J1/J2 < —2.5, the magnetization curve 
exhibits a large jump (of order L'^) on the phase boundary 
between the region of AS^^^ = 1 and that of AS^^^^ — 2, 
whereas for J1/J2 > —2.4 the magnetization curve ap- 
pears to become continuous (i.e., steps are of order L^^). 
Similarly, we observed a large jump in the magnetiza- 
tion between the AS^^^ = 1 and AS^^^ = 3 regions at 
Ji/ J2 = —3.4, but not at other values. However, we find 
that the magnetization curve at J1/J2 = —2.4 also de- 
velops a sharper change at M = M^ with increasing the 
system size L, which turns into almost a discontinuous 
jump at L — 168 [see Fig. [5]Jb)]. This may suggest that 



the transition becomes first order at L — > oo. Since we do 
not have a good scheme of extrapolation to L — + oo for 
incommensurate values of M, it is difficult to determine 
the order of the transition from the numerical calculation 
alone. More elaborated treatments, especially analytical 
ones, would be required for resolving this issue. 

For the parameters calculated, J1/J2 > —3.6, we find 
that at the saturation field h = hg the magnetization 
curve approaches M = 1/2 continuously in accordance 
with the previous studies2.i2i§ (see also Note added). We 
note that the square-root singularity 1/2 — M oc (/ig — 
h)^'"^ is commonly expected for a continuous transition at 
the saturation field, where soft excitations are described 
as free hard-core bosons or free fermions. 



V. VECTOR CHIRAL PHASE 

In this section we take a detailed look at correlation 
functions in the vector chiral phase. We will show that 
this phase corresponds to the low-field regime where the 
magnetization curve has ^S^^^ = 1 steps, and the ground 
state exhibits a long-range order of the longitudinal vec- 
tor chirality, (k["^) = {{si x si+n)') ^ {n = 1,2). We 
first give a brief review of a low-energy field theory for the 
vector chiral phase. We then present numerical DMRG 
results and compare them with the theory. 



A. Bosonization approach for |Ji| <^ J2 



A field theoretical approach to the vector chiral phase 
in the 'h-J2 model was developed by Nersesyan et al. 
in Ref. !29, in which the antiferromagnetic J1-J2 chain 
with easy-plane anisotropy was considered. Kolezhuk 
and Vekua extended this theory to include effects of 
the external magnetic field in Ref. 29. Here we follow 
their approach and apply it to the ferromagnetic J1-J2 
chain. The theory is based on bosonization of the an- 
tiferromagnetic Heisenberg spin chain and perturbative 
renormalization-group (RG) analysis valid for | Ji | <C J2 . 

In the limit | Ji| ^ J2, the model ([1]) can be regarded 
as two antiferromagnetic Heisenberg spin chains which 
are weakly coupled by the ferromagnetic interchain in- 
teraction Ji. Therefore, we apply the standard bosoniza- 
tion technique to the two chains separately, treating the 
interchain coupling Ji as a weak perturbation. The 
low-energy physics of the Heisenberg chains (n — 1,2) 
are described by free bosonic fields {4>n,dn) satisfying 
the equal-time commutation relation [(j)n{x) , dyOn' (y)] = 
i5{x — y)5n,n'- The spin operators s; on the site I — 2j + n 
(j £ Z) in the Hamiltonian ([T]) are expressed in terms of 



the bosonic fields as 

- {-iyasin[2TTMj + V4^0„(a;„)] + ■■■ , (11) 



+ 6'e^v^^"(^") sin[27rMj + V4^(/)„(a;„)] + 



(12) 



where a, b, and b' are nonuniversal constants4^i^ We 
have introduced the continuous space coordinate x, on 
which the bosonic fields depend. On the lattice site I = 
2j -I- n the coordinate x takes the value Xi = j — 1/4 and 
X2 ^ j + 1/4. Equations (fTTj) and (fT2|) allow us to write 
the interchain interaction in terms of the bosonic fields. 
The resulting effective Hamiltonian is given by^ 



^ ^ St/ 



dx 



K, 



d0^ 
dx 



1 



dx 



- gi / dx sin(\/ 87r(/)_ -|- nM) 

52 / dx — — sm(v27rt^_) 
J dx 



with 



51 = JiO^ sin(7rM), 



52 



Ji 

2 



27r6^ 



(13) 



(14) 



Here we have introduced bosonic fields for symmetric (+) 
and antisymmetric (— ) sectors, (/)± = (0i±(/>2)/v2, 6± = 
(6'i±02)/'\/2- In lowest order in Ji the TL-liquid parame- 
ters K± and the renormalized spin velocities v± are given 



v± = f I 1 ± Ji 



K 

■KV 

K 

nv 



(15) 
(16) 



where K and v are respectively the TL-liquid parame- 
ter and the spin velocity of the decoupled antiferromag- 
netic Heisenberg spin chains. The TL-liquid parame- 
ter K is a function of M increasing monotonically from 
K{M = 0) = 1/2 to K{M = 1/2) = Ij^Ld^^ In the 
weak-coupling limit the velocity v is of order J2, except 
near the saturation limit M -^ i, where v ^ and the 
bosonization approach breaks down. 

As we can see from Eqs. ([TT|) and P^ . the gi ((72) 
term in Eq. (fT3|) originates from the longitudinal (trans- 
verse) part of the interchain exchange coupling. These 
coupling constants are renormalized as energy scale is 
decreased in RG transformation. The low-energy physics 
of the effective Hamiltonian p3)) is then determined by 
the strongest of the renormalized coupling constants. In 
case the gi term is most relevant, the 0_ field is pinned 
at a value which minimizes gism{\/8TT(j)^ + ttM). The 



resulting ground state is in the nematic phase, as we 
will discuss in Sec. IVII The vector chiral phase arises 
when the 32 coupling is most relevant and renormalized 
to strong coupling first. 

The scaling dimensions of the 171 and 172 terms are equal 
to 2K- and 1 + {2K-y^, respectively at Ji = 0. It 
is then natural to expect that the 52 term can domi- 
nate over the gi term only in high fields [i.e., for K^ > 
(1 + \/5)/4] when | Ji| ^ ^2-^^ However, as we discussed 
in Sec. IIIII the two-magnon pairing is the strongest in- 
stability at h = hg, which favors the nematic order near 
the saturation field (Fig. [1] and Table [l|. In fact, the 
vector chiral phase is found to be realized in the weak- 
field regime where the bare value of the coupling gi is 
very small {M ^1) and where the classical value of the 
vector chirality is larger {6'^ « '^/'^)- 

In the following discussion let us assume that the renor- 
malized 32 is the largest coupling. In this case we may 
employ the mean-field decoupling scheme introduced by 
Nersesyan et al.^ whose conclusions have been con- 
firmed by numerical studies i21i2^i^i^iH In this scheme 
we assume that both d9+/dx and sin(v27r0_) acquire fi- 
nite expectation values so that the 172 term is minimized. 
We have essentially two choices: 



and 



d0+ 
dx 



/d0+ 
\ dx 



-{TT-2Q), (17a) 



{tt~2Q), (17b) 



where Q is an incommensurate wave number in the trans- 
verse spin correlation (0 < Q < tt/2); see Eq. (pTjl . Note 



that 



and (dO^/dx) have the same sign in the frus- 



trated ferromagnetic chain Ji < 0. The Z2 symmetry is 
spontaneously broken when the ground state selects one 
of the two choices in Eq. p?)) . 

Once the mean-field decoupling is made, the excita- 
tions in the antisymmetric sector ((^_ , 0_) acquire a finite 
energy gap, and the field 0_ , which is dual to the pinned 
field 6-, fluctuates strongly. Therefore the gi term can 
be safely ignored. The symmetric sector (0+,0+) is gov- 
erned by the Gaussian model. 
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dx 
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(18) 



once we redefine the 6+ field, 0+^0+ — {dO+/dx)x, to 
absorb the nonvanishing average {dd+/dx). Hence the 
ground state is a one-component TL liquid. 

We are now ready to calculate correlation functions. 
Most important of these is the ground-state average of 
the vector chirality (O, 



Ji)\ 



-b^{sm{^/2^e^)) ^ Tb^ci, 



/d9. 



(19a) 



= -\l-C2(^-^) = TC2in-2Q), (19b) 



where ci and C2 are positive constants. These nonva- 
nishing averages indicate that the ground state breaks 
a Z2 symmetry and has a vector chiral long-range order. 



Since 



-(i)\ 



and 



,(2)\ 



have the same sign and satisfy Eq. 



([5]), the spin current J„(kj ) flows as depicted in Fig. [31 
The spin current circulates in each triangle in alternating 
fashion, and there is no net spin current flow through the 
whole system. 

Two-point correlation functions are also calculated us- 
ing Eqs. ([n|), nil), ini), and ^. Here we remind the 
reader our convention that the site index I in the orig- 
inal lattice Hamiltonian ^ is equal to 2j + n, where 
the integer j is the site index in each antiferromagnetic 
Heisenberg chain (n = 1, 2) in the two-chain (zigzag lad- 
der) picture. Hence Al = 2Aj — 2A2;. The correlation 
function for the vector chirality is given by 



Kq Ki ) 



- /^(2)\2 



1- 



1 



K+[iTr-2Q)l] 



+ 



(-l)'^K 



cos(2iiMl) 



(20) 



where A^, is a constant (oc b'^gf). In the two-point func- 
tion of the chiral operator kj the uniform 1/P term in 
Eq. ((20I) is replaced by a l/l'^ term^ The transverse and 
longitudinal spin correlation functions are obtained as 



A 



^0^1 . 



{SnS 



\l\l/iK, 



■ cos{Ql) 



0^1 , 



Ar 



K. 



(21) 
(22) 



Lastly the nematic correlation function shows a faster 
decay. 



A' 



[Sq S^ S; Sj_|^-^; 



\l\^/^ + 



cos(2Ql) 



(23) 



In the above equations A and A' are nonuniversal con- 
stants, and we have omitted subleading algebraically- 
decaying terms and exponentially-decaying terms, such 
as a short-ranged incommensurate correlation [ex 
cos(7rM/)] in Eq. (HH). 

We note that the wave number of the transverse spin 
correlation function - the pitch angle - is shifted from 
the commensurate value 7r/2 to Q. The vector chiral 
long-range order and the incommensurate transverse spin 
correlation are the hallmark of the vector chiral phase. 



B. Numerical results 

Here we present our numerical results, which support 
the theory of the preceding subsection. The calculation 
was done for finite open chains with L = 96 and 120 
sites, unless otherwise mentioned. In the following, we 
show mainly the results for L = 120, while we note that 
the results for L = 96 exhibits essentially the same be- 
haviors as those for L = 120. The number of kept DMRG 
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FIG. 10: (Color online) Averaged correlation functions 
for L = 120 spins with J1/J2 = -2.7 and M = 0.1 in 
the vector chiral phase; (a) vector chiral correlation func- 



.(i)^(i)\ 



,(2)^(2) 



(1)^(2)\ 



,(") 



tionS (Kq Kr )av, (Kq Kr )av, and (Kq «r )av, whcrC K- 

(sr X Sr+n)^ [Eq. ([S}] , (b) transverse spin correlation function 






r + l 



(s§s^)av, (c) nematic correlation function 
Open symbols represent the DMRG data. Truncation errors 
are smaller than the size of the symbols. The solid lines and 
solid circles in (b) and (c) are fits to Eqs. (|2H) and (|23|) . re- 
spectively. 



states is up to 350. We have performed typically 10-30 
DMRG sweeps in the calculation and checked the conver- 
gence of the results. Since Eqs. (fT9|) - ((23)) are obtained 
for infinite-length chains, we need to take care of open- 
boundary effects in the DMRG data to make meaningful 
comparison between the theory and the numerics. To 
reduce the boundary effects, we calculate two-point cor- 
relation functions for several pairs of two sites {I, V) with 
fixed distance r = |/ — /'j selecting the two sites being as 
close to the center of the chain as possible. We then take 
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FIG. 11: Averaged, normalized vector chiral correla- 
tions at the distance r — L/2 for L = 120 spin zigzag 



chain. Circles, squares, and triangles represent {kq 



(1)^(1) 



(2J2/lJil)2{«;f 4%). and (2J2/I Ji|){4^'4%), respectively, 
(a) J1/J2 = -2.5, (b) J1/J2 = -2.6, (c) J1/J2 = -2.7, and 
(d) J1/J2 = —2.8. Solid lines are the guide for the eye. Verti- 
cal dashed lines represent the boundaries of the chiral ordered 
phase in low magnetization regime and the chiral disordered 
phase in the high magnetization regime. The shaded region 
corresponds to the magnetization jump at the first-order tran- 
sition. 



their average for the estimate of the correlation. We use 
the notation (••■)av for the averaged correlation func- 
tions below. 

Figure fTOFa) shows a typical r dependence of the av- 
eraged vector chiral correlation functions in the vector 
chiral phase (J1/J2 = -2.7 and M = 0.1). We clearly 
see that the vector chirality is long-range ordered i^ We 



but 



also find that not only (kq Kr )av and (kq k,- ) 
also (kq Kr- )av are positivc, in agreement with Eqs. ^ 
This indicates the ferro-chiral order as drawn in Fig. [3l 

In Fig. [11] we plot averaged vector chiral corre- 
lations at a distance r — L/2 normalized by the 



,(1)^(1) 



(1)^(2) 



and 



coupling ratio, (^ ^K^/a), (2^2/|^i|)(k^ k^/s) 

(2J2/|Ji|)2(4^^4%)' as functions of M and J1/J2. The 
three quantities agree, as expected from Eq. ([5]). This 
figure clearly shows where the vector chiral correlation 
is strong; The vector chiral order exists in the low-field 
regime, but disappears in the high-field regime. The 
shaded regions in Fig. [TT] correspond to the magneti- 
zation jump discussed in Sec. IIVI where the transition 
is clearly first order. We note that near the boundary 
of the vector chiral phase the vector chiral correlations 
suffer boundary effects in open chains and may underes- 
timate the vector chiral order (see the discussion at the 
end of this section). 

We fitted the DMRG data of the transverse spin and 
nematic correlation functions to Eqs. (PTjl and (P5)l . tak- 
ing K^, Q, and the amplitudes A or A' as fitting param- 
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FIG. 12: (Color online) J1/J2 dependence of (a) the expo- 
nent K+ and (b) the incommensurate wave number Q for 
L = 120 spin zigzag chain in the vector chiral phase. Solid 
and open symbols in (a) represent the estimates obtained 
from the fitting of (soS^)av and {sQ'sJs~s7+i)av, respectively. 
The error bars represent the difference of the estimates ob- 
tained from the fitting of the data of different ranges. In (b), 
only the results from (s5s^)av are shown since the estimates 
from {si)Sr)av and {so^s^s^s~_^j)av are identical to each other 
within their error bars. Dashed line in (b) represents the 
classical pitch angle, (ff = arccos(— J1/4J2). 



eters. In the fitting procedure, we used the data of the 
averaged correlation function (sgs^)av [('So''S5^s~s~_^]^)av] 
for L/12 < r < L/2 {L/6 < r < L/2). Figures [TOIb) 
and (c) demonstrate good agreement between numerical 
data and the fits to Eqs. ^ and ^. Note that both 
correlation functions are incommensurate. 

The exponent K^ and the incommensurate wave num- 
ber Q obtained from the fitting are shown in Fig. [T^ 
The TL-liquid parameter K^ is found to be in the range 
0.4 < K+ < 0.8 for the cases we examined numerically, 
and the transverse spin correlation function (sgs^)av is 
the most slowly decaying correlation function except for 
the long-range ordered vector chirality. This suggests 
that a magnetic spiral long-range order in the plane per- 
pendicular to the magnetic field should be realized in real 
three-dimensional materials with additional weak inter- 
chain couplings. The wave number Q, which represents 
the incommensurability of the transverse spin correla- 
tions, shows little dependence on M and decreases as 
J1/J2 decreases. The classical analogue of the wave num- 
ber Q is the pitch angle (j)'^ = arccos(— J1/4J2), which 



shows qualitatively the same feature, but takes a much 
smaller value. We have thus found that the incommen- 
surate wave number is highly renormalized by quantum 
fluctuations towards the commensurate value 7r/2. 

Strictly speaking, the bosonization theory of the previ- 
ous subsection is not directly applicable when the inter- 
chain ferromagnetic coupling is strong, | J1I/J2 ^ 1. The 
consistency between the theory and numerics, as demon- 
strated by the successful fitting, can be understood, once 
we postulate that the vector chiral phase extends to the 
limit J1/J2 — *■ 0, where the bosonization approach is 
valid, and that the low-energy physics in this phase is 
governed by the same effective theory. 

In passing we note that for J1/J2 ^ —2.9 we have not 
found clear evidence for long-range order of vector chi- 
rality; although the vector chiral correlation is strong, it 
seems to decay slowly at long distances in finite-size sys- 
tems L < 120. A possible explanation for this behavior 
would be that the energy gap in the antisymmetric sector 
((/)_, 6'_) is so small that the correlation length becomes 
very large. As a result, the bending-down behavior of the 
vector chiral correlation functions,^" observed near the 
open boundaries (r ^ L) in Fig. fTCTa). penetrates into 
the bulk region and spoils long-range order. For clarify- 
ing the fate of the vector chiral order for J1/J2 ^ —2.9, 
calculations for much larger systems are needed. That is 
left for future studies. 



VI. NEMATIC PHASE AND SPIN DENSITY 
WAVE PHASE 

In this section we discuss in detail the nematic phase 
and the SDW2 phase, where two magnons form a bound 
state with total momentum k = tt. As we saw in Sec. IIVI 
the total magnetization changes in units of AS'^qj = 2 as 
a result of simultaneous fiip of two spins forming a bound 
state. We first review the application of the bosonization 
theory described in Sec. IV Al to these phases, for the sake 
of completeness of our discussion which partly comple- 
ments earlier worksi "^'^'^^ We then present an alternative 
phenomenological theory^ in which a bound magnon pair 
is regarded as a hard-core boson. The theoretical picture 
developed in these discussions is subsequently confirmed 
by numerics. 



A. Bosonization theory revisited 

The nematic phase can be described within the 
bosonization approach for | Ji| ^ J2. When the gi cou- 
pling in the effective Hamiltonian P^ is the most rele- 
vant, the field 0_ is pinned at a value which minimizes 
the gi term. 



.«)^y|Q-M 



(24) 
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In this case the uniform part of the dif feren ce of two 
neighboring spins vanishes, sf — s 



i+i 



0, 

indicating that two spins are boundi^i^ii'^S The dual field 
9- is strongly fluctuating, and the (72 is irrelevant. The 
fields (/>+ and 6+ of the symmetric sector remain gapless 
and constitute a one-component TL liquid P^ . 

The long-distance asymptotic form of correlation func- 
tions can be readily obtained from Eqs. pT|) . (H^), ([M]) . 
and (jlSp . We briefly summarize the results below. 

The longitudinal spin correlation has an incommensu- 
rate oscillatory component, 

(sgsf ) ^M'-^ + jj^ cos [^Z Q - M 

(25) 

where i? is a positive constant (ex a^) and subleading 
terms are omitted here and in the equations below. The 
third term in Eq. (|25p represents incommensurate spin 
density wave correlation. 

The transverse spin correlation {sq s'j^) is short-ranged, 
whose correlation length is the inverse of the gap in the 
{(j>-,9-) sector. Physically, this gap corresponds to the 
binding energy of the two-magnon bound state. 

The composite operator sfsj.^ creating a two-magnon 
bound state represents the nematic order. The nematic 
correlation is alternating and quasi-long-ranged. 
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with B' and B" are positive constants. Here we note that 
two down spins are created and annihilated at neighbor- 
ing sites in Eq. (pS)) . In fact, the algebraic decay with 
the same exponent can be obtained as long as both dis- 
tances between the created down spins and between the 
annihilated ones are odd integers. However, when these 
separations are even, the nematic correlation functions, 
i^o ^tn^7 ^7+2n') {'n^n' <C Oi ^'^^ expected to be weaker, 
as they involve the gapped 9- field in lowest order. This 
is in accordance with the observation we made in Sec. IIIII 
that the wave function of the two-magnon bound state 
with k — IT &i a. saturation field is a linear combination of 
the states in which the distance between two down spins 
is restricted to odd integers. 

Comparing Eqs. (|25p and (I26p . we find that the ne- 
matic correlation (j26p is the most dominant correlation 
if Kj^ > 1. This phase is called nematic phase. On the 
other hand, if K^ < 1, the most dominant correlation 
is the longitudinal spin correlation ((25)l . In this case we 
have the SDW2 phase. 

The vector chiral correlation k^'^' shows a power-law 
decay, 



(2) (2) 



K+P 



(27) 



The same 1/P decay (with a different prefactor) is 
expected for (kq K; ), as the operator product of 
sin(\/27r0_) and the irrelevant g2 term in the Hamilto- 
nian Ti generates the d9+/dx operator. 

In his pioneering paper, ^ Chubukov suggested that the 
nematic phase should have spontaneous dimerization, 
(^f (*f+i ~ ^f-i)) ^ (~1)'- However, in the bosoniza- 
tion theory the dimerization operator is proportional to 
cos(\/87r(/)_ -I- ttA/), whose average vanishes because of 
Eq. ([M)) . We thus conclude that the nematic phase does 
not have a spontaneous dimerization. 



B. Hard-core Bose gas of bound magnons 

As we discussed in Sec. IIIII when —2.7 < J1/J2 < 0, 
the fully polarized state becomes unstable as a result of 
formation of two-magnon bound states at the saturation 
field hs ,L2i^i^ Below hg , bound magnon pairs collectively 
form a TL liquid with nematic correlation as well as in- 
commensurate, longitudinal spin correlation. Here we 
develop a phenomenological theory for the phases which 
emerge as a result of proliferation of p-magnon bound 
states with momentum fc = tt, by assuming that tightly 
bound p-magnons can be treated as a hard-core boson^ 
and ignoring internal structure of the bound states. This 
is expected to be a good approximation as long as the 
density of hard-core bosons is very low, i.e., near the 
saturation field. As we see below, for the p = 2 case, 
this theory is equivalent to the {(f>^,9^) sector of the 
bosonization theory. 

We denote creation and annihilation operators of a 
hard-core boson by b] and 6;. Under the assumption 
that p magnons are tightly bound, we may relate the 
creation operator and density operator of bosons to spin 
operators. 
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(28) 
(29) 



where the (—1)' factor in Eq. (P5)l is introduced because 
the total momentum of bound p magnons is fc = tt; see 
Sec. IIIII In Eq. (pS)) we identify the site index I of the 



boson creation operator b^ with the center-of-mass coor- 



dinate of bound magnons, I = 
we find the density of bosons. 



-(p-l)/2. FromEq. 



-K^" 



(30) 



At small but finite density < p <C 1 the hard-core 
bosons are a TL liquid at low energy. ^^ Its low-energy 
effective theory is again a free field theory, 
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K \dx 



(31) 
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where the bosonic fields (</), 6) play the same role as the 
((/)+, 6*+) fields in the p = 2 case. Here, we take the 
lattice spacing between l-th and (/+l)-th sites to be unity 
and identify I (or I) with x. The TL-liquid parameter 
K depends on interactions that work between bosons in 
addition to the short-range hard-core repulsion. When 
hard-core bosons are free, K ~ I. This should be the 
case in the low-density limit M ^ ^ . In the continuum 
limit the operators bj and bjbi are written as^ 



b] = V-Pe^^<^ Y. 



^2in{-!Tpx + y^4>) 



(32) 



l,\bi^ p+ —=^ + p cos{2TTpx + %/47r</)) H . (33) 

a/tt dx 



as the nematic TL liquid extends from the saturation 
limit (M — > i) to the weak inter-chain coupling regime 
|Ji| < J2; see Fig. m 

To compare the above theoretical results with numeri- 
cal data from DMRG calculation, we need to modify Eqs. 
(IM)l and ([55]) to include finite-size and boundary effects. 
This can be done by calculating the correlation functions 
with Dirichlet boundary conditions on (^{x). Here we can 
borrow results of such calculations from Refs. |4^ and|4y, 
in which correlation functions of the spin-i XXZ model 
in magnetic field are obtained for open spin chains of 
length L, once we notice the mapping of the hard-core 
boson system onto the spin-i XXZ chain [5*;" = (— 1)'6; 

and Sf = bjbi — ^]. In this way we obtain local spin 
polarization (sf) in the J1-J2 spin chain of length L, 



Using these bosonization formulas, it is straightforward 
to calculate correlation functions of hard-core bosons. 
Equations (|^ and (|^ allow us to express these correla- 
tion functions with the original spins si ; we thus obtain 
the longitudinal-spin and p-magnon (multipolar) corre- 
lation functions, {sfsf,) and (s^ • • • s^ ^s^T • • • S;T_,_ j^), 
from the density-density correlation function and the 
propagator of the bosons, respectively. In the thermo- 
dynamic limit L — > cx) we find 
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where 



{st)^^{l-p)-pz{l;q), 



q (-l)'sin(g/) 
z(i\q) = a — ^ — , 
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(36) 

(37) 
(38) 
(39) 



('34') The site dependence of the polarization (s^) comes from 



Friedel oscillations at open boundaries. Such oscillations 
are absent under periodic boundary conditions. The 
characteristic wave vector "2kF" in the Friedel oscilla- 
tions is found from Eqs. ^, ^, and ^ for L > 1 
to be 



where A^, A^, and A^ are positive constants, and 
the parameter 77 in the exponents is related to the 
TL-liquid parameter K hy r] = 2K. Since creating 
less than p magnons costs a finite energy, we expect 
that the transverse-spin correlation functions (sfsf,) and, 
more generally, (s/" • • • s^^^pi^iS^ ■ ■ ■ s'jj+p'-i) with p' < p 
should be short-ranged. 

When p = 2, Eqs. ^^ and ^5^ coincide with Eqs. 
(PS)) and (pS)) by using the relation (|30p and setting 
the exponent rj — K^. That is, the two theoretical 
approaches, the weak-coupling bosonization theory^ for 
|Ji| <C J2 and the phenomenological hard-core boson 
theory^ for i — M <C 1, give a consistent description of 
the nematic and SDW2 phases. This is in fact expected. 



2kF — 2'Kp, 



(40) 



which is determined by the density of the hard-core 
bosons, and is inversely proportional to the number p 
of magnons forming a bound state [see the relation ((30|) ]. 
This result can be easily checked by DMRG calculation. 
The longitudinal spin correlation function in the finite 
J1-J2 spin chain is given by 
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where 
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From Eqs. ([M]) and (|¥l]) . the correlation of longitudinal 
spin fluctuations is obtained as 



'(■'i' 



Sf)(sf,)=p2[^(;^^,.g)_^(;.^)^(^,.^)]^ (44) 



We have two unknown parameters, a and 77, in these for- 
mulas, which can be obtained by fitting numerical data 
to these analytical forms. 

Similarly, the multipolar correlation is obtained as 
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(45) 



which corresponds to the first term in the right-hand side 
of Eq. ([35]) . Here fitting parameters are ^m and 77. 

Before closing this section, we note once again that 
the phenomenological hard-core boson theory is applica- 
ble to any phase which appears as a result of the for- 
mation of p-magnon bound states with p — 2, 3,4, • • • 
and k — t:. We will show in the subsequent sections 
that the phenomenological theory gives a good descrip- 
tion of correlation functions not only in the nematic and 
SDW2 {p = 2) phases, but also in the triatic and SDW3 
(p = 3) phases, and quartic {p = 4) phase, which appear 
for larger |Ji|/J2- Another advantage of the theory is 
that it gives a clear intuitive picture of low-energy exci- 
tations. It is however unable to describe correlations that 
are related to internal structures of bound states [such as 
the feature we discussed below Eq. ([26)) ]. 



C. Numerical results 

We apply the phenomenological hard-core boson the- 
ory to analyze numerical results in this section. We 
use the DMRG method to compute the longitudinal and 
transverse spin correlations, the nematic correlation func- 
tion, and the local spin polarization in finite open chains. 
We fit the correlation functions to Eqs. (P51) - (|15)l with 



p — 2, taking the exponent rj and the coefficients a and 
A,„ as fitting parameters. Since the formulas already in- 
clude effects of open boundaries, here we do not have to 
take spatial average of the correlation functions in the fit- 
ting procedure. Finite-size effects are also properly taken 
into account in these formulas. Indeed, we have observed 
that fitting of numerical results for L = 96 and L = 120 
yields the same good quality of agreement between the 
numerical data and the fits, and gives essentially the same 
estimated values of the fitting parameters. The results 
for L — 120 are shown below. 

Figures [T^ a) and (b) show DMRG results of (sf) and 
{sfsf,)-{sf){sf,) calculated for J1/J2 = -2.0 at M = 0.2 
and 0.4. Shown in the same figures are the fits to Eqs. 
([5S|) and pi)) . respectively, with p — 2. We have used 



and data 
and 



numerical data for 6 < Z < L — 5 to fit (s 
for n < \l - l'\ < L - 10 to fit {sfsf,) - {sf){i^ 
obtained excellent agreement for both. (Note that there 
are only two free parameters rj and a in the fitting.) The 
data of (sf) show Friedel oscillations whose wave length 
is in good agreement with the theoretical prediction with 
p = 2 (without any fitting parameter). This is another 
evidence of the formation of bound magnon pairs. We 
observed this consistency in the whole region of the ne- 
matic and SDW2 phases. The results clearly indicate 
that the low-energy physics in this parameter range is in- 
deed described by the effective theory of hard-core bosons 
of bound magnon pairs. 

For the nematic correlation (s^ s^j^SjTs^T, j^), we fit the 
DMRG result to Eq. (liSjl with p ~ 2, ignoring the ad- 
ditional oscillating component seen in the DMRG data 
which would correspond to the subleading term in Eq. 
(|i5)) . We see in Fig. [T3{c) that the leading power-law 
decaying behavior of (—1)'"' (s/"s^^S;7s;7_|_i) is fitted 
rather well by Eq. ((45l) . Figure [iSjd) shows that the 
transverse-spin correlation function decays exponentially, 
as expected from a finite energy cost for breaking a two- 
magnon bound state. These results also support the 
validity of the effective theory of hard-core bose gas of 
bound magnon pairs. 

Figure[T3]shows the estimate of 77 obtained from the fit- 
ting of (siSi,). As M increases, the exponent rj increases 
across the dashed line 77 = 1. Therefore, the ground state 
undergoes a crossover from the low-field SDW2 phase, 
where the longitudinal spin correlation function is dom- 
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FIG. 13: (Color online) Correlation functions for L — 120 
spin zigzag chain with J1/J2 = —2.0, M = 0.2 and 0.4, and 
their fits to the theory for the nematic and SDW2 phases: (a) 
Friedel oscillations in the local spin polarization (sf), (b) lon- 
gitudinal spin fluctuation (sfsf/) — (sf)(sf/), (c) nematic corre- 
lation function (s/^Sj'liSjTSjT , j^)- The open symbols represent 
the DMRG data. Truncation errors are smaller than the size 
of the symbols. In (b) and (c), the data for I = L/2— [r/2] and 
I' = L/2 + [{r + l)/2] are plotted as a function of r = |Z - Z'|. 
The results of the fitting are shown by solid symbols in (a) 
and (b) and by dashed curves in (c). The data for M — 0.4 
are multiplied by a factor 2 in (a) and shifted by 0.01 in (b). 
(d) Absolute values of the averaged transverse-spin correla- 
tion function {sgs^)av. 
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FIG. 14: (Color online) M dependence of the exponent i] in 
the nematic phase (77 > 1) and SDW2 phase (r; < 1). The 
estimates are obtained from the fitting of (sfsi/). The error 
bars represent the difference of the estimates obtained from 
the fitting using the data of different ranges. 



inant, to the high-field nematic phase, where the ne- 
matic correlation dominates »^ We have found that 77 esti- 
mated from the other correlators, (sf), (sfsf,) — (sf){sf,), 



'I I 



^i^v 



and (s; Sj^j^Sj, S;,^j^), are consistent with Fig. [141 As 

Af — ^ i , ?7 increases towards ry = 2, in agreement with 
the theoretical prediction that the hard-core bosons be- 
come free in the dilute limit. 

VII. INCOMMENSURATE NEMATIC PHASE 

As we discussed in Sec. |llTl for -2.720 < J1/J2 < 
—2.669 the fully-polarized state becomes unstable at the 
saturation field as a result of formation of two-magnon 
bound states with an incommensurate momentum.^ Be- 
low the saturation field these bound states are expected 
to form a TL liquid with incommensurate nematic cor- 
relation. Such an incommensurate nematic phase was 
predicted by Chubukovji who dubbed this phase the chi- 
ral biaxial spin nematic, as he considered it to have long- 
range vector chiral order. However, our numerical results 
indicate that the vector chirality is not long-ranged in the 
incommensurate nematic phase. 

Figure [12] shows our DMRG results of spatially aver- 
aged correlation functions for J1/J2 = —2.7 and M — 
0.4. We see in Fig. [TSTa) and (b) that the correlation of 
the longitudinal spin fluctuations [{s^sD — {sf^) (sr)]av and 
the nematic correlation function (S(j"sf s,7s^_|_]^)av decay 
slowly (presumably algebraically) with an incommensu- 
rate modulation. On the other hand, we find in Fig.[T5rc) 
that the transverse-spin correlation (soS^)av decays expo- 
nentially, indicating the existence of a finite binding en- 
ergy of the two-magnon bound pairs. Finally, Fig. llSf d) 
shows that the correlation functions of vector chirality 
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FIG. 15: (Color online) Averaged correlation functions for 
L = 120 spin zigzag chain with J1/J2 = —2.7 and M = 0.4 in 
the incommensurate nematic phase; (a) longitudinal spin fluc- 
tuation [(sqS^) — (so)(Sr)]av, (b) nematic correlation function 
{sq's]^s~s~, j^)av, (c) absolute values of the transverse spin 
correlation function {sQS^)a,v, and (d) vector chiral correlation 



,(i)^(i)\ 



J2)(2) 



(1)^(2)\ 



functions (kq Kr )av, («o '^r )av, and (kq kJ. )av Trunca- 
tion errors are smaller than the size of the symbols. Insets 
in (a), (b), and (d) show the absolute values of the data in a 
log-log scale. 
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FIG. 16: (Color online) Squared modulus of the Fourier 
transform, |s^(fc)| , which is an even function of k and where 
s"(fc) = (l/\/i)Eie''='((sf) - M), for J1/J2 = -2.7 and 
M — 0.4 in the incommensurate nematic phase. The solid, 
dashed, and dotted lines represents the data for L = 160, 120, 
and 80, respectively. Inset: the local spin polarization (sf) 
for J1/J2 = -2.7, M = 0.4, and L = 160. 



K; have at most quasi-long-range order with incommen- 
surate oscillations. 

In Fig. [ini we show the local spin polarization (sf) 
and its Fourier transform for L ~ 80, 120, and 160. The 
Fourier transform exhibits three peaks. This is in con- 
trast to the cases of the nematic phase in Sec. VI and 
the triatic and quartic phases (see Sec. IVIli]) . where the 
polarization (sf) is described by Eq. (|36p with a single 
wave number q. The positions of peaks in the Fourier 
transform are almost independent of L, suggesting that 
the incommensurability should not be due to finite-size 
nor open-boundary effects. 

We have observed qualitatively the same behaviors 
of the correlation functions and spin polarization for 
J1/J2 = -2.7 and M > M^ w 0.35. From these 
results, we conclude that the incommensurate nematic 
phase (without chiral long-range order) exists in the nar- 
row region of the phase diagram; see Fig. [TJ 

Unfortunately, we are not aware of an effective theory 
which can give consistent description of these numerical 
results. However, in the spirit of the hard-core boson 
theory in Sec. VI, we may try to treat the two-magnon 
bound states as hard-core bosons: 
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(46) 

(47) 



where b[ ; and b!^ ; are creation operators of hard-core 
bosons with momentum k — tt + S and tt- J, respectively. 
The long-range order of vector chirality would follow if 
the average of the boson number difference, b\bi — 63^21 
became nonvanishing spontaneously. Figure I15f d) indi- 
cates that this is not the case, and implies that the low- 
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energy properties in this phase are determined by two 
bosonic modes (6{6;^ -^ '>' 
two-flavor TL Hquid. 



&2^2 8'Il'i ^1^1 



62&2); forming a 



VIII. TRIATIC AND QUARTIC PHASES 

In this section we consider the triatic, SDW3, and quar- 
tic phases, in which the total magnetization changes by 
^'5'tot — 3 and 4. We again apply the hard-core boson 
theory of Sec. IVIBI with p = 3 and 4. 

Figure flTl shows our DMRG results for the triatic and 
SDW3 phases. As a typical example, we chose the cou- 
pling ratio J1/J2 — —3.0, and magnetization per spin 
M = 0.2 and 0.4. We fit the local spin polarization 
(sf) and the longitudinal spin fluctuation correlation 
(sfsf') - (sf>(sf'> to Eqs. dSS]) and ^ with p = 3, 
respectively, taking 77 and a as fitting parameters. We 
find that these correlators are fitted quite well by the 
formulas. The fitting of the triatic correlation function 
{s^ s'l'^^s^_^2^i' ^'i'+i^i'+2) to Eq. ([^5]) with p = 3 also 
works well, within the approximation that the sublead- 
ing oscillating terms are ignored. The transverse spin 
and two-magnon (nematic) correlation functions, (sfsf,) 



and 



{S, Sj^-^S, 



^I'+l 



I, decay exponentially, in accordance 



with the theoretical prediction. All these observations 
demonstrate the validity of the bosonic effective theory 
for the triatic/SDWa phase. 

Figure fTSl shows the exponent 77 obtained from the fit- 
ting of (sfsf,) for J1/J2 > -3.0 and {sf) for J1/J2 < 
—3.2. Although the estimates have rather large error 
bars, there is a clear tendency that 77 increases from 77 < 1 
to 77 > 1 as M increases. The estimates of ry obtained 
from the other correlators, including the triatic correla- 
tion (s/'s^]^Sjt|_2S/7s;7_|_]^S;7_|_2) for J1/J2 > —3.0, exhibit 
essentially the same feature. The result indicates that 
the ground state undergoes a crossover from the low- 
field SDW3 phase with the dominant longitudinal-spin 
correlation to the high-field triatic phase where the three- 
magnon (triatic) correlation is dominant. The behavior 
of 77 at large M is also consistent with the theoretical 
prediction that 77 — > 2 as Af ^ ^ . 

For the quartic phase, we show the local spin polariza- 
tion {sf) calculated at J1/J2 = —3.6 and M = 0.2 and 
0.4. (The numerical data of other correlation functions 
are not available for J1/J2 < —3 because of slow conver- 
gence of DMRG calculationi^) We clearly see in Fig. [19] 
that the numerical data of (sf ) are fitted well by Eq. ((5^ 
with p = A. This gives strong support for the presence of 
the quartic phase for these parameters from the follow- 
ing reason. As we emphasize below Eq. (l40l) . the period 
of the Friedel oscillations in (sf) is directly related to 
the density p of hard-core bosons and, in particular, the 
number p of magnons forming a bound state. Indeed, if 
we compare Figs.fTSlal.fTWa). and [T9l for the same mag- 
netization M, we find that the number of nodes in (sf) 
change as ex 1/p for the nematic/SDW2, triatic/SDWs, 
and quartic phases. Therefore, the successful fitting of 




FIG. 17: (Color online) Correlation functions for L — 120 
spin zigzag chain with J1/J2 — —3.0, M — 0.2,0.4, and their 
fits to the theory for the triatic phase: (a) Friedel oscilla- 
tions in the local spin polarization (sf), (b) longitudinal spin 
fiuctuation (sfsf,) — {sf){s^,), (c) triatic correlation function 
('s/^s^-^s,"t_2Sj7Sj7 I jSjT I 2). The open symbols represent the 
DMRG data. Truncation errors are smaller than the size of 
the symbols. In (b) and (c), the data for I — L/2 — [r/2] and 
I' = L/2 + [{r + l)/2] are plotted as a function oir ^\l-l'\. 
The results of the fitting are shown by solid symbols in (a) and 
(b) and by dashed curves in (c). The data for M = 0.4 are 
multiplied by a factor 2 in (a) and shifted by 0.01 in (b). (d) 
Absolute values of the averaged transverse-spin and nematic 
correlation functions, (sgs^)av and (S(J's]'"s~s~^j)av. 
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FIG. 18: (Color online) M dependence of the exponent r\ in 
the triatic phase (77 > 1) and SDW3 phase [j] < 1). The esti- 
mates are obtained from the fitting of {s\s\,) for ,h/ J2 > —3.0 
and {si) for J1/J2 = —3.2, —3.4. The error bars represent the 
difference of the estimates obtained from the fitting using the 
data of different ranges. 



(sf ) to Eq. (pS)) with a certain p can be considered as 
a strong evidence for the formation of p-magnon bound 
states. 

One may naturally expect that for the quartic phase a 
spin-density- wave (SDW4) regime with dominant longi- 
tudinal spin correlation should also appear at low mag- 
netic field. Indeed, the exponent 77 obtained from the 
fitting of (sf) for J1/J2 = —3.6 (not shown here) ex- 
hibits a tendency that 77 changes from 77 < 1 to 7/ > 1 
with increasing M . Meanwhile, the estimates of -q have 
large error bars, which prevent us from determining accu- 
rately the crossover point between the quartic and SDW4 
regions, unfortunately. We therefore tentatively call the 
region of /S.S^^^ — 4 the quartic phase, keeping it in mind 
that the low-field part of the phase probably includes the 
SDW4 regime. 



IX. CONCLUDING REMARKS 

We have determined the magnetic phase diagram of 
the spin-1/2 J1-J2 zigzag spin chain with ferromagnetic 
Ji and competing antiferromagnetic J2 interactions un- 
der magnetic field. Asymptotic behaviors of correlation 
functions have been derived for various phases, using 
bosonization approach for weak Ji and the effective the- 
ory for hard-core bosons of bound multi-magnons. By 
fitting numerical data of the correlation functions ob- 
tained by the DMRG method to the analytic forms, 
we have successfully identified the vector chiral phase, 
nematic/SDW2 phases, triatic/SDWs phases, and quar- 
tic phase. 

At low magnetic field, we have found the vector chiral 
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FIG. 19: (Color online) Friedel oscillations in the local spin 
polarization (sf ) for L = 120 spin zigzag chain with J1/J2 = 
—3.6 and M — 0.2, 0.4, and their fits to the theory for the 
quartic phase. The open symbols represent the DMRG data 
and solid symbols show the results of the fitting. Truncation 
errors of the DMRG data are smaller than the size of the 
symbols. The data for M — 0.4 are multiplied by a factor 2. 



phase, marked by long-range vector chiral order and alge- 
braically decaying incommensurate transverse-spin corre- 
lations. The vector chiral state is the quantum counter- 
part of the classical helical state. In the classical Ji- J2 
model the helical state appears as the ground state in 
the whole magnetization region, whereas in the spin-1/2 
model the vector chiral phase appears only in the low- 
magnetization regime (but not at M = 0). For larger 
magnetization, the chiral state is destroyed by the for- 
mation of quantum magnon bound states, and turns into 
the spin density wave states. 

At higher magnetic field, magnons form stable bound 
states, caused by ferromagnetic attractive interactions. 
The stabilization of magnon bound states is a general fea- 
ture of frustrated ferromagnets. Spin multipolar orders 
induced by the bound-state formation were discovered 
recently in spin-1/2 models on the square lattice;^ the 
triangular lattice;^ and the two-leg ladder latticei^ In 
all of these models, magnon bound states become stable 
when the ferromagnetic state is destroyed by compet- 
ing antiferromagnetic and/or ring-exchange interactions. 
The unique feature of the present zigzag chain model 
is that the number of magnons forming a bound state 
increases consecutively with approaching the ferromag- 
netic phase boundary J1/J2 = —4. This unique feature 
might relate to the fact that the zigzag spin chain with 
J1/J2 = —4 has highly degenerate ground states;'^^''^^ 

We have found various phases that can be well de- 
scribed by the effective theory for hard-core bose gas of 
bound multi-magnons. Near the saturation field, there 
appear various multipolar TL liquid phases, such as the 
nematic (quadrupolar) , triatic (octupolar), and quartic 
(hexadecapolar) phases, which are characterized by the 
condensation of bound multi-magnons. With lowering 
the magnetic field from the saturation field, the increase 
of bound-magnon density enhances the effect of repul- 
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sion between bound magnons. This leads to a crossover 
from a high-field region of the multipolar TL liquids to a 
low-field region of spin density wave states, where density 
waves of bound multi-magnons dominate. 

Lastly we note that our phase diagram appears to be 
qualitatively consistent with the magnetic properties ob- 
served in the spin-1/2 chain cuprate LiCuV04, whose 
exchange couplings are estimated as J1/J2 ~ —0.4 with 
ferromagnetic Jif^i Recent experiments revealed that, 
when a magnetic field is applied, the low-temperature 
phase undergoes two successive transitions with changing 
fLeldji2,iii Besides an usual spin-flop transition at hd ~ 
2.5T, which is presumably due to spin anisotropy effect, 
another magnetic phase transition occurs at hc2 ~ 7.5T. 
Below hc2, (or more precisely, in hd < h < hc2), the 
low-temperature phase has spiral spin structure, having 
incommensurate spin order in the plane perpendicular 
to the applied field, as well as ferroelectlicityJ^ii^ii^ In 
the field above hc2, the system has a modulated mag- 
netic order parallel to the field while the perpendicular 
spin components are disorderedJ^ In comparison with 
the magnetic phase diagram of the J1-J2 model, it is nat- 
ural to identify the magnetic transition at /i = /ic2 with 
the transition between the vector chiral phase and the 
SDW2 phase in our model ([T]). This means that the ex- 
perimentally observed spin modulated state in high field 
can be characterized by a density wave order of bound 
magnon pairs. Furthermore, in the light of our phase 
diagram, we predict that the nematic phase, which has 
not yet been observed experimentally, should appear in 
higher magnetic field. Further experimental studies on 
high magnetization phases would be interesting. On the 
theoretical side, it is important to include effects of inter- 
chain couplings, further interactions, and spin anisotropy 
to make more quantitative comparison between our the- 
oretical results and experiments. In LiCuV04, the inter- 
chain coupling indeed induces three-dimensional order at 
very low temperature T < Tn ~ 2.3K. We also note that 
the critical field hc2 is about 0.2 of the saturation field 
hs w 41T and considerably larger than the value esti- 
mated from the one-dimensional J1-J2 model. 

Note added: Since the submission of this paper, a 
preprint by Sudan et al.^^ has appeared, in which a phase 
diagram very similar to ours is obtained independently. 
They have found a direct metamagnetic transition from 
the vector chiral phase to the ferromagnetic phase for 
large IJ1I/J2 where bound states oi p > 5 magnons are 
formed. In Sec. IIVI we discussed magnetization curves 
for J1/J2 > —3.6 where only bound magnons of p < 4 
participate in the magnetization process. 
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APPENDIX: BLOCK THEOREM FOR SPIN 
CURRENT 

In this appendix we prove that there is no net spin 
current flow in the ground state (even in the vector chiral 
phase). This is a variant of Bloch's theorem that there is 
no net current flow in the ground state without external 
field.^^^ 

We begin with defining the spin current. Since the 
Hamiltonian ([T]) conserves the z component of the total 
spin J2i^h the s^ current is a well-defined quantity. The 
equation of motion for sf reads 



dt'' 



-i[sin] 



-JAn 



.(1) 



.«)-j2(.p^-.[!y,(A.i) 



from which we deduce the s^ current flowing from the site 
I to the site l + nis given by J; = Jn^^i {n = 1,2). It 
also follows from Eq. (|A.1[) that 

m<l 



(2) 



,.(2) 



(A.2) 



implying that the total s^ current flowing through the 
system is Jtot = Jik^^^ + 2J2k'-^\ where we have sup- 
pressed the site index I. 

To prove our statement that the net current Jtot always 
vanishes in the ground state, we consider the J1-J2 spin 
chain of finite length L with the periodic boundary con- 
dition and add to its Hamiltonian Ti. a weak symmetry 
breaking term 



I 



(A.3) 



with the coupling y > 0. This term allows us to select a 
unique ground state with a broken Z2 symmetry. Let us 
assume that the unique ground state \g)L,y has a nonva- 
nishing expectation value of a linear combination of the 
vector chiralities. 



JlK 



(1) 



2J2K1 



L,y 



> 



3y|Ji 



(A.4) 



where (• • • )L^y is the average in the ground state \g)L 



We now introduce the twist operator 



V 



Uo 



exp 



1=1 



(A.5) 
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with the twist angle = 2tt/L. We then take Ue\g)L,y as 
a trial state, which has a smaller net spin current than 
\g)L,y, and compare its energy with the energy of the as- 
sumed ground state \g)L,y We find the energy difference 

AE = {Ulin + ny)Ue)L,y " (H + Hy)L,y 

= Ji(cos6i - 1) ^(sf sf+1 + s^s^^^ - ynf^)L,y 
I 

+ J2(COS20 - 1) ^(sf sr+2 + sf sf+2>L,y 
I 

- Jisin0^(Kp)+y(sfsf+i +sfsf+i))i,j, 



J2sin20^(KpOi,y 



(A.6) 



Assuming the translation invariancc of the ground state, 
we reduce Eq. (|A.6p to 

AE = -27r(JiKp^ + 2.hnf'^)L,y 

- 27ryJi (sf.r+i + ^hl+i)L,y + ^(^"') (A-7) 

for L 3> 1. Using the inequality —Ji{sfs'f^-^ + 
sfsf+i)L,a < fl-^ili we conclude from Eqs. (IA.4P and 



(jA.yp that AE < 0. This is in contradiction with the 
assumption of \g)L.y being the ground state. This means 
that our assumption (jA.4[l is not valid, and instead we 
have 



(JlK[')+2J2Kp^)L,y<ay, 



(A., 



where a is a positive constant (a = 3| Ji|/4). 

In the same way, starting from the assumption that the 
ground state \g)L,y has a negative net current, (JiK; + 

2J2K1 )L,y < —jy\Ji\, and using the twisted trial state 
Ue\g)L.y with the angle 9 — — 27r/L, we can again show 
that the trial state has a lower energy than the assumed 
ground state, and thereby we have (Jik[ +2J2k\ iL.y > 
—ay. We thus obtain 



.(1 



(2) 



\{Jii^r+2J2nY')L,y\<ay 



(A.9) 



We now take the limit L — > 00 and then y — > 0, yielding 
Ji(Kp)) + 2J2(Kp^)=0. (A.IO) 

It is straightforward to generalize this identity to the case 
when the ground state breaks translation symmetry as 
well as to other spin Hamiltonians. 
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